c
c -- file name gx2phs.htm 030100 C SUBROUTINE GXDISP is called from Group 13, only if the C PATCH name starts with KEDI. The PATCH type should be CELL, C and COeff = GRND1, VAL = 0 for both KE and EP for the Mostafa C & Mongia model; COeff = GRND2, VAL=0.0 for the Chen and Wood C model; and COeff = FIXFLU, VAL = GRND3 for the bubble-induced C turbulence production model ( for more details see PHENC entry C 'TURBULENCE MODELLING FOR TWO-PHASE FLOWS'). c SUBROUTINE GXDISP INCLUDE 'farray' INCLUDE 'grdear' INCLUDE 'grdloc' INCLUDE 'satgrd' INCLUDE 'satear' common/genl/dbfil1(14),debgtz,dbfil2(45) COMMON/NAMFN/NAMFUN,NAMSUB CHARACTER*6 NAMFUN,NAMSUB EQUIVALENCE (JPRODB,EASP3) LOGICAL EQZ logical dbfil1,debgtz,dbfil2 C IF(IXL.LT.0) RETURN NAMSUB='GXDISP' IF(INDVAR.EQ.KE.OR.INDVAR.EQ.EP) THEN IF(ISC.EQ.2.OR.ISC.EQ.3) THEN C.... Put large-scale turbulent motion timescale Te into EASP2 CONST=0.35 IF(IGR.EQ.3) CONST=0.165 CALL FN15(EASP2,KE,EP,0.0,CONST) IF(STORE(LBNAME('TE '))) CALL FN0(LBNAME('TE '),EASP2) C.... Put characteristic particle response time Tp into EASP3 CALL FN56(EASP3,VOL,DEN2,INTFRC,1.0) CALL FN26(EASP3,R2) IF(STORE(LBNAME('TP '))) CALL FN0(LBNAME('TP '),EASP3) C.... Put equivalent of 2.RHO2.R1.R2./Tp into CO CALL FN21(CO,R1,INTFRC,0.0,2.0) C.... CO = GRND1 - Mostafa & Mongia IF(ISC.EQ.2) THEN C.... Put Te+Tp into EASP8 CALL FN10(EASP8,EASP2,EASP3,0.0,1.0,1.0) C.... Now get Tp/(Te+Tp) and put in EASP3 CALL FN27(EASP3,EASP8) C.... Now get CO = CO * (F(EASP3+I) = 2.0*R1*R2*FIP*(Tp/(Te+Tp)) C (I from 1 to NY*NX). CALL FN26(CO,EASP3) C.... CO = GRND2 - Chen & Wood ELSE IF(INDVAR.EQ.KE) THEN L0CO=L0F(CO) L0TE=L0F(EASP2) L0TP=L0F(EASP3) DO 20 IX=IXF,IXL DO 20 IY=IYF,IYL I=IY+NY*(IX-1) F(L0CO+I)=F(L0CO+I)*(1.0-EXP(-0.0825*F(L0TP+I) 1 /(F(L0TE+I)+TINY))) 20 CONTINUE ENDIF ENDIF C.... VAL = GRND3 - Bubble-induced turbulence production ELSEIF(ISC.EQ.15) THEN IF(INDVAR.EQ.KE) THEN IF(EQZ(EL1A)) EL1A=0.05 CALL FNVSLP(1,L0F(JPRODB),CFIPA) cc call prn('pr 1',jprodb) CALL FN55(VAL,INTFRC,JPRODB,JPRODB,EL1A) CALL FN0(JPRODB,VAL) IF(STORE(LBNAME('PRKB'))) CALL FN0(LBNAME('PRKB'),JPRODB) cc call prn('pr 2',jprodb) ELSEIF(INDVAR.EQ.EP) THEN CALL FN56(VAL,JPRODB,EP,KE,C1E) cc call prn('epvl',val) ENDIF CALL FN26(VAL,R1) ENDIF ENDIF NAMSUB='gxdisp' END c C**** SUBROUTINE GXIPST is called from GREX3, group 13, and is entered C only when the patch name begins with character 'IPST'. C SUBROUTINE GXIPST(INTFRC,CINT,NPATCH) INCLUDE 'farray' INCLUDE 'grdloc' INCLUDE 'satgrd' COMMON /IGE/IXF,IXL,IYF,IYL,IREG,NZSTEP,IGR,ISC,IRUN,IZSTEP,ITHYD, 1 ISWEEP,ISTEP,INDVAR,VAL,CO,NDIREC,WALDIS,PATGEO,IGES20(6) COMMON/IDATA/NX,NY,IFL2(118) COMMON/GENI/NXNY,IGFIL(59) CHARACTER*(*) NPATCH INTEGER VAL,CO,WALDIS,PATGEO LOGICAL FLUID1,CNV,IPT,NEZ,GTZ COMMON /NAMFN/NAMFUN,NAMSUB CHARACTER*6 NAMFUN,NAMSUB NAMSUB = 'GXIPST' C IF(CNV(INDVAR).AND.IPT(INDVAR) .AND. CINT.GT.0.0) THEN C.... Obtain zero-location indices L0IP = L0F(INTFRC) L0VAL = L0F(VAL) L0VAR = L0F(INDVAR) L0VOLD = L0F(OLD(INDVAR)) IF(FLUID1(INDVAR)) THEN C.... Variable belongs to first phase INDPRT = INDVAR + 1 L0U = L0F(U1) ELSE C.... Variable belongs to second phase INDPRT = INDVAR - 1 L0U = L0F(U2) ENDIF L0PRT = L0F(INDPRT) L0POLD = L0F(OLD(INDPRT)) FACTOR = 0.25*CINT IF(NPATCH(5:6).EQ.'FE') THEN C.... Fully explicit DO 10 I = 1,NXNY F(L0VAL+I) = 0.0 IF(NEZ(F(L0U+I))) THEN INEXT = I + NY IF(GTZ(F(L0U+I))) INEXT = I - NY F(L0VAL+I) = FACTOR*F(L0IP+I)* 1 (-2.0* (F(L0PRT+I)-F(L0VAR+I))+F(L0POLD+I)- 1 F(L0VOLD+I)+F(L0POLD+INEXT)-F(L0VOLD+INEXT)) ENDIF 10 CONTINUE ELSEIF(NPATCH(5:6).EQ.'CN') THEN C.... Crank-Nicholson FACTOR = 0.5*FACTOR DO 20 I = 1,NXNY F(L0VAL+I) = 0.0 IF(NEZ(F(L0U+I))) THEN INEXT = I + NY IF(GTZ(F(L0U+I))) INEXT = I - NY F(L0VAL+I) = FACTOR*F(L0IP+I)* 1 (-3.0* (F(L0PRT+I)-F(L0VAR+I))+ 1 F(L0PRT+INEXT)-F(L0VAR+INEXT)+F(L0POLD+I)- 1 F(L0VOLD+I)+F(L0POLD+INEXT)-F(L0VOLD+INEXT)) ENDIF 20 CONTINUE ELSE DO 30 I = 1,NXNY F(L0VAL+I) = 0.0 IF(NEZ(F(L0U+I))) THEN INEXT = I + NY IF(GTZ(F(L0U+I))) INEXT = I - NY F(L0VAL+I) = FACTOR*F(L0IP+I)* 1 (-F(L0PRT+I)+F(L0VAR+I)+F(L0PRT+INEXT)-F(L0VAR+INEXT)) ENDIF 30 CONTINUE ENDIF ENDIF NAMSUB = 'gxipst' END c